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present an implementation and discuss first test applications to electron and ion transfers in 
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Keywords: Density functional theory, Car-Parrinello molecular dynamics, external potentials 



I. INTRODUCTION 



Electrical fields may induce charge transport or also 
electrochemical reactions. To simulate such condensed- 
phase phenomena, it is desirable to have such an exter- 
nal field implemented in a code which allows to perform 
Car-Parrinello dynamics pQ using periodic boundary con- 
ditions [2]. The theoretical treatment of external fields 
within Kohn-Sham theory [3 , 4 is in principle straightfor- 
ward. The practical implementation in a code employing 
periodic boundary conditions, however, leads to the ques- 
tion how potential discontinuities of the potential energy 
at the borders of the unit cell should be treated. Sev- 
eral approaches to tackle this problem have been devised 
in the past years [5TTLQ]. In the present paper we use a 
somewhat different approach with the aim to devise an 
implementation in the CPMD code which is both concep- 
tually simple and simple to use in practical applications, 
even for time-dependent external fields. We achieve this 
by representing an external field by sine functions where 
the problem of discontinuities does not arise. 
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The energy is calculated with the Kohn-Sham expres- 
sion [3 0]: 
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We add the time-dependent external potential to the 
Kohn-Sham energy expression: 



II. THEORETICAL BACKGROUND 

Car-Parrinello molecular dynamics is based on an ex- 
tended Lagrangian pQ: 
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This results in the following equations of motion: 
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We have implemented into the Car-Parrinello molecu- 
lar dynamics code [2] a sinusoidal potential field which 
may change temporally and spatially in three dimensions. 
For one spatial dimension: 



^sineO, t) = (a(t) + A) • sin ( 27r— + b(t) 

' or. 



a(t), b(t) may be linear or sinusoidal functions. With 
this flexible implementation, a broad range of phenom- 
ena can be simulated. However, it must be considered 
that the wave length cannot exceed the size of the unit 
cell as the spatial part 2izx/l x ensures that the periodic 
boundary conditions are fulfilled. 
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FIG. 1. Arrangement of the gold-dibenzenedithiol-gold junc- 
tion in the supercell. (grey: C, white: H, yellow: S, gold: Au) 



III. APPLICATIONS 
A. Electron Motion: Time-Dependent Amplitude 

Metal- molecule-metal junctions are used to measure 
the electron transport across single molecules. Many 
experimental and theoretical studies aimed at under- 
standing the mechanism of the electron transport (see 
[TT] and literature cited therein). It was found that 
the resulting current depends on many factors which 
are difficult to define experimentally. In particular the 
precise arrangement of a single molecule between two 
surfaces or tips is not easy to control in experiment. 
Since most of the molecules under investigation are 
isolators and the distance between the surfaces or tips is 
in the order of a nanometer, it is not clear from experi- 
ment if the effective current is significantly influenced by 
the linking molecules. Theoretical investigations were 
usually performed using perturbation theory approaches 
only. Dynamics calculations with an explicit field give 
a more complete picture of the evolution of a complex 
system and may be used in the future to study all 
these questions in detail. In the present study we 
want to illustrate in first generic test calculations the 
influence of a three-dimensional electrical field on a 
metal- molecule-metal junction. 

We studied a system consisting of a dibenzenedithiol 
molecule covalently linked to two gold clusters. The 
electric field is applied along the junction (Figure 1). 

A field with the following parameters was used: 
A = 0.000025 a.u. 

/ 27T 7r\ 
a(t)= 0.000025- sin * - - j a.u. 

b(t) = 

The amplitude varies sinusoidal with time between 
and 0.00005 a.u. (0.03 kcal/mol). The transfer of one 
electron along the molecular junction is observed. 

The snapshots in Figure 2 show the motion of a 
localized (Wannier) orbital. The total spin density is 
essentially zero during the full process. The external field 
causes a shift of the total charge rather than a build-up 
of spin density. This shift occurs in a continuous motion. 
The total electron density is just slightly shifted to one 



40.7 fs, 0.000019 a.u. 




48.5 fs, 0.000025 a.u. 




53.3 fs, 0.000029 a.u. 




FIG. 2. Snapshots of the motion of a localized orbital during 
the application of an external field with increasing amplitude 
to a gold-dibenzenedithiol-gold junction. 

side. 

Figure 3 shows the course of the energies and tem- 
perature during the simulation run. The total energy is 
raised by about 40 kcal/mol due to the application of 
the external field. 

The electrons follow essentially instantaneously even 
if the fictitious mass of the electrons is high in the Car- 
Parrinello molecular dynamics. The process is essentially 
adiabatic and is determined by the change of the ampli- 
tude of the external field. These first results indicate that 
the junction acts as a nearly perfect insulator and could 
just as well be omitted to obtain the same effect. 

B. Ion Migration: Time-Dependent Phase 

/3-Eucryptite (LiAlSiCU) is known as a one- 
dimensional Li + ionic conductor [T2] , The conductivity 
is along the crystallographic c-axis of the quartz-like 
structure. In the CPMD simulation one crystal defect 
is introduced in the supercell (Figure 4). A Li + ion is 
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Amplitude 0.000000 0.000005 0.000017 0.000033 0.000045 0.000050 a.u. 
Time 0.0 19.4 38.7 58.1 77.4 96.8 fs 



FIG. 3. Total energy of the Car- Parr inello Lagrangian (blue), 
energy of the classical Hamiltonian (green), Kohn-Sham en- 
ergy (black), temperature (red) during the application of an 
external field to a gold-dibenzenedithiol-gold junction. 





FIG. 4. Unit cell of /3-eucryptite. Within the unit cell the 
lithium ions may migrate along four different channels. Dark 
blue: Si, blue: Al, red: O, gold: Li, yellow: Li channel with 
vacancy. 

removed and a Si 4+ ion is substituted by an Al 3+ ion to 
obtain a neutral system. The potential field is applied 
along the c-axis. The phase of the potential field varies 
linearly in time. 

A field with the following parameters was used: 

A = 0.0001 a.u. 

a(t) = 



The field which has an amplitude of 0.0001 a.u. (0.06 
kcal/mol) is moving with a velocity of 15000 m/s, 
respectively. This corresponds to a phase shift of 2tt 
within 0.29 ps. The series of snapshots in (Figure 5) 
shows the migration of five Li + -ions to new lattice sites 
on a time scale of a few hundred femtoseconds. The 
motion is started by a lithium ion near the left border of 



FIG. 5. Snapshots from a CPMD simulation showing the 
motion of lithium ions (yellow) in /3-eucryptite. The channel 
containing a vacancy is shown in light yellow. Dark blue: Si, 
blue: Al, red: O. The arrows mark the migrating lithium ions. 
The migration can also be followed from the trajectories to 
the left. Five lithium ions in the left half of the simulation 
cell shown migrate to the next lattice site generating a new 
vacancy at the left border of the simulation cell. 

the simulation cell. The neighbouring lithium ions follow 
till the original vacancy is filled and a new vacancy near 
the left border of the simulation cell is generated. On a 
longer time scale this vacancy would be filled by lithium 
ions from the right border (respecting periodic boundary 
conditions). 

In Figure 6 the temporal evolution of the energies and 
temperature is shown. The total energy shows a slight os- 
cillation in the beginning when the field is turned on. The 
kinetic energy of the electrons (difference between the to- 
tal energy and the energy of the classical Hamiltonian) 
stays small. The temperature (corresponding to the dif- 
ference between the energy of the classical Hamiltonian 
and the Kohn-Sham energy) shows a slight increase by 
about 30 K. Obviously the ion migrations are possible 
with a very small amplitude of the electric field and a 
moderate take-up of energy. 

IV. CONCLUSIONS 

We have added a sinusoidal field to Car-Parrinello 
molecular dynamics. The field can change with phase 
and amplitude. Like this it is possible to model a vast 
diversity of experimental situations in which an electri- 
cal field matters. The setup can easily be extended to 
more complex periodic fields by using more than one si- 
nus function. In first test calculations we have illustrated 
the application to electron and ion transfers. 
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FIG. 6. Total energy of the Car- Parr inello Lagrangian (blue), 
energy of the classical Hamiltonian (green), Kohn-Sham en- 
ergy (black), temperature (red) during the application of a 
temporally moving external field to /3-eucryptite. 



plane wave code [2 j. For all calculations the BLYP 
exchange-correlation functional [T3l [14] was used in its 
unrestricted formulation (LSD). Troullier-Martins pseu- 
dopotentials were employed for describing the core elec- 
trons [15j. The pseudopotential cutoff was set to 70.0 Ry- 
dberg. The fictitious electron mass was set to 400 a.u. 
and a time step of 4 a.u. (0.097 fs) was used. After the 
equilibration, the temperature was not controlled. 
For the simulations of the metal- molecule- metal junc- 
tions a periodically repeated simulation cell with a size 
of 15 x 15 x 60 A 3 was used. The systems were initially 
equilibrated at a temperature of 150 K in order to reduce 
molecular vibrations and to focus on the electronic mo- 
tion. 

For the simulations of the eucryptite crystal periodic 
boundary conditions were applied with a cell size of 10.5 
x 9.1 x 44.8 A 3 . The systems were equilibrated at 300 K. 



METHODS 



The simulations were performed using Car-Parrinello 
molecular dynamics pQ as implemented in the CPMD 
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